Pathlength distributions of atmospheric neutrinos 



T.K. Gaisser & Todor Stanev 
Bartol Research Institute, University of Delaware, Newark, DE 19716 

We present the distribution of the production heights of atmsopheric neutrinos as a function 
of zenith angle and neutrino energy. The distributions can be used as the input for evaluation of 
neutrino propagation under various hypotheses for neutrino flavor oscillations. 



I. INTRODUCTION 

■ Initial results from SuperKamiokande appear to confirm indications from 1MB, Q Kamiokande (|] and Soudan jj] 
t-H ' of an excess of v e relative to in the atmospheric neutrinos. One possible interpretation is that neutrino flavor 
^q. oscillations play a role. In a two-flavor mixing scheme, for example, the probability that a neutrino of flavor i and 
^ ' energy Ei retains its identity after propagating a distance L in vacuum is Q 
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where dm 2 is the difference in mass squared of the two neutrino mass eigenstates and 9 is the mixing angle. Therefore, 
to evaluate the manifestation of the mixing in a detector that measures to some degree the direction and energy of 
neutrino-induced events, one needs to know the distributions of production heights of the neutrinos as a function of 
energy and zenith angle. More complicated mixing schemes [[| and effects of propagation in matter still require 
this basic information about the points of origin of the neutrinos. 

Information about origin of the neutrinos is implicit in any calculation of neutrino fluxes. Here we extract the 
relevant information from the simulation of Rcf. |^| , which has been compared to several other calculations in Ref . [|| . 

The paper is organized in three sections. First we review the simulation we are using to calculate production 
of neutrinos in the atmosphere. Next we present the basic results of the calculation. We discuss simple analytic 
approximations which offer insight into the systematics of the results and compare them to simulation results for 
Oh' zenith angles from the vertical to horizontal. Finally, we provide some paramctrizations, based on the analytic 
q approximations, that may be useful for practical application of the results. 

C3 ■ II- SIMULATION 



The simulation was performed in the spirit of earlier calculations of the atmospheric neutrino flux |^,|To|. The 
simulation code is one dimensional. In this approximation, all secondaries are assumed to move in the direction of 
the primary particles (except for a small fraction of low energy secondaries with angles larger than 90° to the beam, 
which are discarded). The validity of this approximation has been checked in Refs. JTT|^2[ . 

The primary cosmic ray flux and its composition is the parametrization used previously in the calculation of Agrawal 
et al. H which in the multi-GeV range falls in between the measurements of Refs. [ p^JT^ ]. Incident cosmic-ray nuclei 
are treated in the superposition approximation Ga|, with cascades generated separately for protons and neutrons in 
order to insure the correct ratios of neutrinos and antineutrinos. The fraction of neutrons is derived from the fractions 
of nuclei heavier than hydrogen in the primary flux. 

We consider three ranges of neutrino energies that correspond approximately to the three major types of experi- 
mental events in a detector the size of SuperKamiokande: contained events; partially contained neutrino interactions 
and stopping neutrino induced muons; and througoing muons. The energy ranges are presented in two different ways: 

• 0.3 < E v < 2 GeV; 2 < E v < 20 GeV; E v > 20 GeV and 

• E > 1, 10 and 100 GeV. 

The integral form is more closely related to simple analytic approximations that we use as the basis of parametrizations 
of the results of Monte Carlo simulations. 

Our results are obtained with the geomagnetic cutoffs for Kamioka and for the epoch of solar minimum, which is 
applicable to measurements performed currently (~ 1994 — 99). Because of the high geomagnetic cutoffs at Kamioka, 
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it is not necessary to account precisely for the phase of the solar cycle. To illustrate the potential influence of 
geomagnetic effects at other locations we also tabulate some results for the much higher geomagnetic latitude of the 
SNO experiment. 

We have not included prompt neutrino production through charm decay because it is totally negligible in the 
considered energy ranges [^^]. All neutrinos are generated either in pion and kaon decays or in muon decays. The 
production heights are stored separately for neutrinos from n/K and from muon decays. The muon decay procedure 
accounts for the muon enegy loss during propagation in the atmosphere. Technically the muon lifetime is sampled 
in the muon rest frame and then the muon is propagated in the atmosphere with time dilation proportional to its 
decreasing energy. Thus mouns decay on the average sooner that they would have if one (incorrectly) sampled from 
a decay distribution using their Lorentz factor at production. 



III. RESULTS 



Before presenting the results for neutrinos, we show a comparison between measurements of GeV muons at different 
altitudes in the atmosphere and our calculation made with the same Monte Carlo code [jl9). This type of balloon 
measurement provides the most direct test of the validity of the cascade model and of the treatment of the muon 
propagation in the atmosphere because the muons and neutrinos have a common origin. The comparison shown in 
Fig. 1 is with data of the MASS experiment JH| as discussed in Ref. fl9|| . 

Fig. 2 shows the height of production of neutrinos of energy above 1 GeV for cos(#) = 0.75. The graph gives dN v /dh 
(cm _2 s _1 sr _1 km _1 ), where h is the slant distance from the neutrino production point to sea level. Contributions 
from muon decay and from tt/ K decay are shown separately for v e + v e and for + T>^ . The overall flux of v e + z7 e 
from ir/K decay is much lower because it reflects primarily the contribution of K\ decays, which is very low in this 
energy range. The curves for electron and for muon neutrinos from muon decay are nearly equal. They extend to 
lower altitudes with a slope that depends on the average energy of the parent muons. For higher energy neutrinos 
this slope is signicantly flatter as a consequence of the higher parent muon energy and correspondingly longer muon 
decay length. For E v > 20 GeV, most parent muons reach the ground (except in nearly horizontal direction) and stop 
before decaying. As a consequence, the height distributions for neutrinos from muon decay deep in the atmosphere 
are nearly flat. 



A. Height distribution for neutrinos from tv/K decay 
1. Analytic approximation 

It is instructive to look at a simple approximation for the height of production of neutrinos from decay of pions. 
In the approximation of an exponential atmosphere with scale height Hq and the approximation of Feynman scaling 
for the production cross sections of pions in interactions of hadrons with nuclei of the atmosphere, a straightforward 



solution of the equations for propagation of hadrons through the atmosphere fl2(T gives 21 



dX An 7(7 + 1) 

for the integral flux of neutrinos in the energy range E v <C where reinteraction of pions in the atmosphere can be 
neglected. There is a similar expression for neutrinos from decay of kaons proportional to Bk X Znk- The meaning 
and approximate values of the quantities in these equations are given in Table Q. 



TABLE I. Values of the parameters used in Eq. 3.1 that correspond to a power law primary cosmic ray spectrum and to 
an exponential atmosphere. 7 and K are the spectral index and the coefficient of the differential cosmic ray energy spectrum, 
dN/dE = KE~ < -~ t+1 \ 7v (r_ff) is (m^/m^) 2 ((m^/mx) 2 )- Zn-tt (Zmk) is the spectrum weighted moment for pion (kaon) 
production by nucleons (Zn-x = J dxx 1 dN/dx). Ajv and Ajv are the attenuation and interaction lengths for nucleons. Bk 
is the branching ratio for K — > /1 decay. Xq and ho are the total vertical thickness (in g/cm 2 ) and the scaleheight for an 
exponential atmosphere in km. 



7 


K 


Ajv 




Znw 


Bk x Znk 






X 


ho 




cm- 2 S - 1 sr- 1 (GeV) 7 


g/cm 2 


g/cm 2 










g/cm 2 


km 


1.70 


1.8 


86 


120 


0.08 


0.0075 


0.5731 


0.0458 


1030 


6.4 



2 



In Eq. 3.1 X is the slant depth in the atmosphere at which the pion is produced and decays. We now convert this 
into distance I from the detector in the approximation of an exponential atmosphere in which 



X 
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(3.2) 



For 9 < 70°, curvature of the earth can be neglected and cos 9 in Eq. 3.2 is cosine of the zenith angle to a good 
approximation. The effective values of cos 9 for larger angles are given below. 

The corresponding approximate expression for the distribution of production distances is 
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where X is to be evaluated as a function of £ from Eq. [3.2j . Assuming a primary cosmic ray nucleon flux with the 
normalization given in Table | (and including the small contribution from decay of kaons) the normalization factor is 
AXo/ho ^ 0.020. 

Taking parameters from Table | gives the most probable distance of production as 

k ° ' X ° (3.4) 



which is 



In- 

cos 9 An cos 9 
15 km for vertical neutrinos from decay of pions. 



2. Monte Carlo results 



Fig. 3 shows the distance distribution for neutrinos from tt/K decay [E v > 1 GeV) for cos(#) = 1.00, 0.75, 0.50, 
0.25, 0.15 and 0.05. Here 9 is the zenith angle of the neutrino trajectory at the surface of the Earth. For large zenith 
angles, the curvature of the Earth is significant, and it is necessary to use effective values of cos e ff(9) that represent 
the convolution of the locations of neutrino production with the local zenit h a ngle as it decreases moving upward 
along the trajectory. We treat cos e ff(9) as a free parameter in fitting Eqs. 3J: and 3J3 to the Monte Carlo results. 
The values are included in Table ||. 

TABLE II. Comparison of analytic and Monte Carlo values of the effective value of cos 9 . Column 1 shows the cosine of 
the zenith angle 0. Column 2 shows the most probable production height for neut rinos from tt/K decay from the Monte Carlo 
calculation. Column 3 gives the most probable height of production from Eq. 3.4 with cos e //# from column 4. Columns 4 & 5 
give the cos e // 9 values that fit best the calculated height of production distributions for neutrinos from tt/K and muon decay 
with ho = 6.50 km. Column 6 gives the normalization coefficient C M needed to fit the distribution for neutrinos from muon 
decay. 
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1.00 
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1.00 
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1.00 

0.75 

0.50 

0.26 
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0.087 



0.69 
0.71 
0.77 
0.83 
1.00 
1.86 



Up to cos 9 — 0.25 the agreement between the Monte Carlo calculation and the analytic estimate is quite good. Note 
that for nearly horizontal neutrinos the height distribution from the Monte Carlo calculation is artificially narrow and 
irregular. The atmospheric model used does not treat exactly the atmospheric densities at vertical depths of less than 
few g/cm 2 . This intruduces a sharp cutoff in the height distribution for strongly inclined showers and also decreases 
the width of the height distribution. 

The height distribution of v e from tt/K decay has a similar shape with much lower normalization, because only K\ 
have decay mode with v e 's (K® 3 ). 
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B. Height distribution for neutrinos from muon decay 



1. Analytic approximation 

To estimate the height of production for neutrinos from decay of muons is more complicated because of the com- 
petition between decay and energy loss for muons in the multi-GeV energy range. One starts from the distribution 
of production points for muons, which is similar to Eq. |3.1| with different coefficients. The resulting approximate 
expression pj[ for the distribution of production distances (differential in the energy of the parent muons as well as 
the slant height of production) is 
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where t is the muon lifetime, 
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At high altitude muo n energy loss (a(X — Y)) can be neglected and the expression 3.5 is proportional to slant depth 
X given by Eq. 3.2. This expression gives a good account of the high-altitude exponential falloff of the neutrinos 
from muon decay. 

An approximation that is adequate for fitting the distribution for all distances (integrated over neutrino energy) is 
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where C M is an overall normalization factor used to fit the Monte Carlo results (see Table |J) • 



2. Monte Carlo results 



Fig. 4 shows the height of production distributions for muon neutrinos of energy above 1 GeV from muon decay. 



The lines are calculated accordin g to Eq. 3.6 with values of cos e // 9 as given in Table |I. To obtain the fits shown in 
Fig. 4 the approximations of Eq. |3.6| have also been renormalized as indicated in Table ij. 

At high altitude the height of production distributions have the same shape as the ones from neutrinos from ir/K 
decay, shifted to lower altitudes by one muon decay length (6.24 km for 1 GeV muons). At lower altitude the shapes 
are quite different. The production height for neutrinos from muon decay extend to much lower altitude because of 
the slow attenuation of the parent muon flux, an effect which becomes more pronounced as the energy increases. 

It is interesting to observe that at high zenith angles the yield of neutrinos from (daughter) muon decay exceeds 
the yield of neutrinos from the decay of the parent pions and kaons. The reason is that muon neutrinos from muon 
decay in flight have a spectrum extending almost to x = 1 (where x = E^/En), while the neutrinos from 7r decay 
can only reach E™ ax = E n X (1 — r„) = 0.428 E„. The corresponding Z-factors Z vfxv and Z vv are 0.133 and 0.087 
respectively, including the effect of muon polarization in pion decay. The result is that for large zenith angles, when 
almost all muons decay the yield from muon decay becomes slightly larger than that from n/K decay (~ 9/7 for 
cos(9 = 0.05). 

Fig. 5 compares the distributions of distance to production for from muon decay with E v above 1, 10, and 100 
GeV. At high neutrino energy the muon decay length becomes comperable or larger than the total dimension of the 
atmosphere. The neutrino height of production then becomes constant deep in the atmosphere. 

The height distribution for v e from muon decay is analogous to that of . The only difference is the slightly lower 
normalization, which reflects the ratio Z llVlx jZ ilVe = 0.133/0.129 = 1.03 (including muon polarization). 
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C. Height distribution in three energy bins 



In Tabic III we show the average height of production and the contributions of n/K and muon decays for neutrinos 
in the three energy bins (0.3 < E v < 2 GeV; 2 < E v < 20 GeV; 20 GeV) which roughly correspond to contained 
neutrino events, semicontained events and stopping neutrino induced muons and throughgoing neutrino induced 



muons. For each angle and neutrino flavor Table III first gives the average height of production (slant depth) in km 
and the width of the height of production distribution. Then it gives the contribution (in %) of tt/K decay f m and 
the corresponding (h m ) and ah m , then the same quantities (/ M , (h^), ah^ for neutrinos from muon decay. 

The calculation was done with the geomagnetic cutoffs of Kamioka, except for the three lines (cos 9 — 1.00, 0.75 and 
0.50, for the lowest energy bin) that are also calculated for the high geomagnetic latitude of SNO (^2|. The numbers 
for high geomagnetic latitude are slightly higher (2 - 10 %) for both neutrino sources because of the contribution of 
low energy protons. This difference becomes negligible at higher angles. 

TABLE III. Production height (slant distance, km) of neutrinos for six values of cos 8 and three neutrino energy ranges. 
The calculation is for the geomagnetic location of Kamioka with three lines for the lowest energy range calculated for Sudbury, 
Canada. 
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There are two obvious trends in the numbers in Table III. The height of production for neutrinos from tt/K decay 
grows slightly with the neutrino energy because higher energy mesons preferentially decay (rather than interact) in 
the tenuous atmosphere at high altitude. Neutrinos from muon decay, on the other hand, are generated at lower 
altitude at high energy because of the increasing muon decay length. This second feature is much stronger because 
of the proportionality of decay length and muon energy. 

The average heights of production also reflect the relative yields of the two neutrino sources. For low energy v e [y e ) : 
for example, the contribution of K® 3 is small, so the average height of production is dominated by muon decay. At 
higher energy the relative contribution of K® 3 grows, especially at directions close to the vertical, and (h) becomes 
intermediate between those of the two processes with correspondingly larger width. 

Generally the contribution from muon decay increases significantly with the zenith angle since even 20 GeV muons 
easily decay in cascades developing in nearly horizontal direction. 
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IV. CONCLUSIONS 



We have calculated the distribution of pathlengths of atmospheric neutrinos for a range of angles and energies 
relevant for current searches for neutrino oscillations with atmospheric neutrinos. Accounting correctly for the path- 
lenght will be particularly important for neutrinos near the horizontal direction where the pathlength through the 
atmosphere of neutrinos from above the horizon is of the same order of magnitude as the pathlength through the Earth 
of neutrinos from below the horizon. We have also given simple approximations that may be useful in interpolating 
the tables and adapting the results for different energy ranges and directions. 

The influence of the geomagnetic effects on the calculated height of neutrino production distributions is not very 
strong. The difference in the average production heights for neutrinos detected at Kamioka and SNO is of order 
several per cent in directions relatively close to the vertical. This difference diminishes with angle and becomes totally 
negligible for upward going neutrinos, where the geomagnetic cutoffs becomes approximately equal, being averaged 
over the geomagnetic fields of the opposite hemisphere. 

The agreement of our calculation with the measured muon fluxes above 1 GeV/c as a function of the atmospheric 
depth serves as a check on the validity of the results presented above. 
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FIG. 1. Comparison of the calculated flux of negative muons above 1 GeV/c as a function of the atmospheric depth to the 
measurements of the MASS experiment. ^j|Jl9| 
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FIG. 2. Height of production distribution for neutrinos of energy above 1 GeV at cos 6 = 0.75. The contributions of n/K 
and muon decays for the two neutrino flavors are clearly visible. Dots show the height of production for electron neutrinos 
from n/K decay, short dashes: u e + from muon decay. The heavy dash line is the sum of the two. Dash dot: + z? M from 
tt/K decay, dash-dash: + £/ M from muon decay. The heavy solid line is the sum of these two. 



FIG. 3. Height of productio n for muon neutrinos above 1 GeV from tt/K decay as a function of angle with different symbols. 
Lines show the best fits of Eq. 3.3 with the parameters given in Table O. 



FIG. 4. Height o f pr oduction for muon neutrinos above 1 GeV from muon decay as a function of angle (dots). Lines show 
the best fits of Eq. 3.6 with the parameters given in Table |n| 



FIG. 5. Height of pr oduction for muon neutrinos above 1, 10, and 100 GeV from muon decay at cos# = 0.25. Lines show 
the best fits of Eq. 3.6 and with normalization factors of 0.83, 0.50, and 0.33 for 1, 10, and 100 GeV respectively. 
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Gaisser & Stanev, "Height of production ", Fig. 2 
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Gaisser & Stanev, "Height of production ", Fig. 3 
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Gaisser & Stanev, "Height of production ", Fig. 4 
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Gaisser & Stanev, "Height of production ", Fig. 5 



